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FORCE  FREE  MAGNETIC  FIELDS 
AND  SOLAR  ACTIVITY 

I.  INTRODUCTION 

The  key  phenomenon  of  solar  activity  is  the  solar  flare,  and  it 
is  generally  agreed  that  the  energy  released  in  a solar  flare  is 
derived  from  magnetic  energy  ''Sweet,  1969).  The  structure  of  magnetic- 
field  configurations  in  the  solar  atmosphere  is  therefore  essential 
for  understanding  solar  activity  in  general  and  solar  flares  in 
particular.  One  meaningful  way  to  approach  the  problem  of  flare 


prediction  is  to  develop  methods  for  determining  magnetic-field 
configurations  in  active  regions  and  to  find  procedures  for  predicting 
the  evolution  of  these  configurations.  Since  the  pressure  of  the 
coronal  plasma  is  less  than  the  magnetic  stress  in  active  regions, 
the  magnetic  field  in  solar  active  regions  typically  comprises  one 
of  two  configurations.  One  of  these  is  the  "force-free"  state 
(Schmidt,  1963),  and  the  other  is  the  "current  sheet".  It  is  widely 
believed  that  the  current  sheet  plays  a key  role  in  the  flare  process 
since  the  "Petschek  mechanism"  (Petschek,  196^  constitutes  a rapid 


mode  for  the  release  of  free  magnetic  energy.  Our  studies  of  solar 
flares  ''Barnes  and  Sturrock,  1972)  indicate  that  the  transformation 
of  a force-free  configuration  into  a current-sheet  configuration  may 
be  a key  step  in  the  overall  flare  process. 
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the  free  magnetic  energy  is  the  excess  energy  of  a force-free  field  j . 
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configuration.  In  the  Spicer  model,  rapid  energy  release  becomes  j ae rcia Hrmn n • ^ j 

possible  through  a resistive-type  instability,  comparable  with  the  - — 
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"tearing-mode  instability"  Forth,  Killeen  and  Rosenbluth,  106? ' which 
leads  to  rapid  energy  release  from  current  sheets. 

For  the  above  reasons,  it  is  essential  that  one  be  able  to 
compute  force-free  magnetic-field  structures  if  one  is  to  make  detailed 
analysis  of  magnetic- field  changes  which  lead  to  solar  flares.  This 
is  true  either  for  a short-term  goal  or  a long-term  goal.  The  short- 
term goal  is  simply  that  of  developing  magnetic-field  models  which 
are  simple  enough  for  detailed  study,  but  which  bear  some  relation  to 
magnetic-field  patterns  likely  to  arise  in  the  sun's  atmosphere.  The 
long-term  goal  is  that  of  computing  magne t ic- f ie Id  structures  in  the 
sun's  atmosphere  based  on  study  of  the  observed  magnetic-field  con- 
figuration of  the  photosphere  and  the  photospheric  velocity  field. 


The  latter  involves  the  calculation  of  force-free  magnetic-field 
conf  igurat’’ ons  in  three  dimensions  which  is  beyond  the  capability 
of  the  computer  now  available  to  us.  For  this  reason,  we  have  limited 
our  scope  to  the  study  of  model  force-free  magnetic-field  configurations 
in  two  dimensions. 


II.  SPECIFIC  COALS 


Our  first  goal  was  that  of  developing  computational  procedure 
for  the  calculation  of  force-free  magne t ic- f ie Id  configurations.  The 


relevant  equation  is 


j X B = 


7 x B = 4rrj 
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so  that  equation  (2.1)  may  be  rewritten  as 
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B X (V  X B)  = 0 


(2.3) 


It  follows  from  this  euqation  that 


V X B = kB 


(2.4) 


where  n(x)  is  a scalar  function  of  position.  However,  on  operating  on 
both  sides  of  this  equation  with  the  gradient  operator  and  noting  that 


B = 0 


(2.5) 


we  find  that 


B . V > 


0 


(2.6) 


so  that  m is  a constant  along  any  field  line.  In  general,  however, 
k will  vary  from  one  field  line  to  the  next. 

If  one  considers  the  restricted  class  of  force-free  magnetic -fie Id 
configurations  for  which  k = const.,  the  governing  equation  (2.4)  is  a 
linear  equation.  This  is  a great  simplification  in  computing  field 
configurations  and  it  is  one  which  has  been  used  extensively  in  the 
investigation  of  force-free  fields.  However,  there  is  no  a priori 
reason  to  expect  that  a real  magnetic-field  conf fguration  occurring  in 
the  sun's  atmosphere  will  have  x = constant.  For  instance,  if  we 
consider  the  magnetic-field  configuration  produced  by  a localized 
current  distribution,  we  may  expect  that  it  will  vary  with  distance 
L asymptotically  as  a power  law: 


B * L 
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This  means  that 


(V  X B)  « L‘n-1  (2.8) 

so  that 

| k 1 “ L_ 1 . (2.9) 

It  can,  in  fact,  be  shown  that  a magnetic-field  configuration  constrained 
to  have  k = constant  over  an  infinite  volume  will  have  infinite  energy 
fAnzer,  1976). 

The  procedure  which  we  have  developed  at  Stanford  for  the  study 
of  force-free  magnetic-field  configurations,  without  the  restriction 
k = constant,  is  based  on  the  representation 

B = Vat  x vg  (2.10) 

which  follows  from  equation  (2.5).  One  readily  verifies  that 

B . Va  = 0 , B . 7B  = 0 , (2.11) 

so  that  both  nt  and  B are  constant  along  any  field  line.  Hence  we  see 
also  that 


(7  X Bl  . 79  = 0 , (7  X B)  . 7f  = 0 , (2.12) 

which  may  be  rewritten  as 

7r»  . [7  X (79  X7g)]  = 0 (2.13) 

7B  . [V  X (Va  X VP)]  = 0 . (2.1U) 
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These  may  be  expressed  in  diadic  notation  as  follows 


(V0  )2  - Vpvp ; Wty  - (Vfy.Vp)V' p + VaVp:VVp  = 0 (2.15) 

- (Va.7p)7 2a  + vave:V7a  + (Vn)2v  p - Vo7a; vvp  = o . (2.16) 

These  equations  are  in  fact  suitable  for  numerical  calculation  by 

the  relaxation  procedure.  This  is  equivalent  to  determining  whether, 

— ♦ - 

for  a given  assumed  field  configuration,  j X B vanishes  at  a particular 
point  and,  if  it  does  not,  allowing  the  field  line  to  adjust  itself 
to  make  j X B = 0 locally. 

The  principal  aim  of  this  research  program  has  been  to  study  models 
of  force-free  magnetic-field  configurations  using  the  representation 
(2.1?)  and  the  relaxation  procedure  for  numerical  computation. 

III.  RESULTS 

In  order  to  make  the  models  sufficiently  simple  to  calculate, 
it  is  essential  to  introduce  some  kind  of  symmetry.  It  is  possible 
to  develop  field  models  using  translational  symmetry  in  terms  of 
a rectangular  coordinate  system,  as  was  done  some  years  ago  (Sturrock 
and  Woodbury,  I967 . However,  such  field  patterns  can  never  erupt  to 
form  current  sheets  since  an  open  field  pattern  based  on  a line  dipole 
has  infinite  energy  per  unit  length.  It  is  therefore  essential  to  go 
to  a more  realistic  geometry  such  as  cylindrical  geometry  or  spherical 
geometry. 

When  force-free  field  patterns  of  cylindrical  symmetry  were  first 
investigated  Barnes  and  Sturrock,  1972 'l,  we  found  that  the  total 
energy  in  the  magnetic  field  configuration  increases  progressively 
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with  angle  of  shear,  and  that,  for  shear  angles  above  some  critical 
value  (about  tt),  the  energy  in  the  magnetic  field  exceeds  the  energy 
of  the  corresponding  open-field  configuration.  This  then  led  to  the 
important  possibility  that  a force-free  magnetic-field  configuration 
might,  if  sufficiently  stressed,  become  unstable  against  an  "eruptive 
instability"  which  develops  the  field  pattern  into  an  open  configuration 
involving  a current  sheet.  This  current  sheet  could  then  provide 
further  energy  release  by  field-line  reconnection  (Petschek,  I96U), 
which  is  the  probable  mechanism  for  energy  release  in  solar  flares. 

One  of  the  main  aims  of  this  research  program  has  been  to  further 
investigate  model  force-free  magnetic-field  patterns,  extending  the 
geometry  to  the  more  realistic  spherical  geometry.  It  was  important 
to  determine  whether  the  energy  property  just  mentioned  is  common  to 
all  force-free  configurations  or  whether  it  applies  only  to  some 
special  class  of  such  configurations. 

We  faced  a number  of  problems  in  pursuing  this  aim.  One  real 
handicap  was  the  limited  capability  of  the  Stanford  computer  which 
made  it  necessary  to  work  with  a much  coarser  mesh  than  is  desirable. 
Another  difficulty  was  that  of  obtaining  convergence  of  the  relaxation 
procedure.  It  is  essential  to  introduce  a "relaxation  parameter"  \ 
defined  by 

b'0+1  = (1  - \)  e'n)  + x z(nl  (3.1 

est. 

where  B 0 is  the  field  pattern  obtained  by  the  n'th  iteration  and 
B is  the  field  pattern  indicated  by  "full"  relaxation,  adjusting 

each  mesh  point  in  turn  or  in  blocks)  to  satisfy  the  equations  (2.15) 


# ••  ♦«*  -Vs 


and  (2.16).  As  the  field  became  progressively  more  stressed,  it  was 
necessary  to  use  a smaller  and  smaller  value  of  X.  to  insure  convergence. 
Unfortunately,  this  means  that  the  convergence  became  slower  and  slower. 

Since  there  were  no  simple  analytical  models,  based  on  spherical 
geometry,  against  which  to  test  our  code,  a minor  error  persisted  in 
the  code  for  some  time. 

An  example  of  our  calculation  of  force-free  field  configurations, 
based  on  a spherical  boundary,  is  shown  in  Figures  3.1  - 3*3.  It 
is  assumed  that  the  field  pattern  has  cylindrical  symmetry  above  the 
axis  of  the  spherical  coordinate  system  and  this  leads  to  the  repre- 
sentation (2.10)  for  which  a = a(r,9)  and 

P = CD  - vfr.e)  . (3.2) 

Contours  of  constant  a give  the  field  lines  themselves,  since  a is  a 
measure  of  magnetic  flux  enclosed  in  a circle  with  constant  values  of 
r and  9.  Contours  of  constant  \ give  the  shear  angle  of  the  field 
line.  In  fact,  it  was  convenient  to  use  log  r rather  than  r as  the 
radial  variable.  This  leads  to  a mesh  which  is  finer  in  the  high-field 
regions  and  coarser  in  the  weak-field  regions.  It  was  assumed  that  or 
and  v are  prescribed  along  the  surface  of  the  sphere. 

An  additional  difficulty  in  our  calculations  concerned  the  choice 
of  outer  boundary.  (in  principal,  we  would  have  wished  that  we  could 
dispense  with  an  outer  boundary  and  extend  the  field  pattern  to  infinity.) 
It  was,  however,  convenient  to  introduce  an  outer  boundary  which  was 
perfectly  conducting.  Field  lines  can  therefore  not  penetrate  the 
limits  of  the  computational  mesh.  if  the  outer  radius  of  the  system 


■+>  4 ' * # 
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could  be  chosen  large  enough,  one  could  hope  that  the  conducting  boundary 
would  appear  to  the  field  to  be  lying  at  infinity.  Then  the  exact  form 
of  boundary  conditions  chosen  for  the  outer  spherical  boundary  would  not 
be  important. 

We  found,  however,  that  as  the  field  shear  was  increased,  a point 
was  reached  at  which  the  fields  seemed  to  billow  out  rapidly  (as  a function 
of  rate  of  increase  of  shear,  not  time)  and  to  press  against  the  outer 
boundary.  This  readily  occurred  even  though  the  outer  boundary  was 
placed  as  far  out  as  700  surface  radii.  In  comparing  the  energy  stored 
in  the  fields,  under  various  conditions  of  shear,  with  the  energy  stored 
in  a completely  open  field  configuration  (current  free,  except  for  a 
current  sheet  region),  it  was  observed  that,  as  the  shear  ^fes  increased, 
increasing  the  stored  energy  in  the  system,  the  lines  billowed  out  before 
the  stored  energy  reached  the  equivalent  open  field  energy.  This  implies 
that  there  is  a smooth  transition  between  closed  and  open  field  config- 
urations in  the  spherical  model,  with  no  build-up  of  excess  stored  energy 
as  shear  increases.  Figures  3*1  ~ 3*3  show  some  typical  field  config- 
urations for  surface  shear  angles  of  2.6,  tt  and  2ti  radians,  respectively. 
The  boundaries  are  at  0 = 0 and  tt/2  and  R/Ra  = 1 and  180.  The  boundary 
at  R/Rq  = lBO  is  taken  to  be  open  (approximated  by  allowing  field  lines 
to  pass  through  the  boundary  with  the  condition  that  they  be  straight 
there).  The  equatorial  boundary  (0  = tt/2)  is  taken  to  be  conducting, 
as  if  a mirror  image  existed  in  the  southern  hemisphere.  One  can  see 
in  Figure  3»1  maximum  shear  = 2.6  radians)  that  the  field  lines  bulge 
outward  to  a fairly  large  radius,  but  do  not  seem  to  sense  the  presence 
of  the  outer  boundary  the  field  configuration  for  this  shear  angle 

8 


* ' v 


is  observed  to  be  the  same  for  both  open  and  closed  outer  boundaries). 
Figure  3*2  shows  the  fields  for  a shear  of  tt  radians.  Here,  the  fields 
extend  beyond  the  boundary.  Then,  in  Figure  3-3  (shear  = 2tt  radians), 
the  fields  extend  outward  even  further.  The  total  stored  energy 
increases  only  slowly  with  increasing  shear,  appearing  to  approach 
asymptotically  the  open  field  value. 

Similar  calculations  on  systems  which  include  the  southern 
hemisphere  as  well  show  only  qualitative  differences,  mainly  that  the 
fields  extend  downwards  into  the  southern  hemisphere.  The  outward 
bulging  characteristic  remains  the  same  however. 

The  important  comparison,  for  our  purposes,  is  shown  in  Figure 
in  which  the  energy  of  the  magnetic  field  is  plotted  as  a function  of 
shear  angle  y^  (the  maximum  value  of  y).  One  of  the  most  important 
results  of  our  program  was  the  discovery  that,  with  spherical  geometry, 
the  behavior  is  quite  different  from  that  which  we  previously  found 
with  cylindrical  geometry.  Whereas,  with  cylindrical  geometry,  the 
total  energy  rose  above  the  value  for  the  corresponding  open-field 
configuration,  this  was  not  the  case  for  spherical  geometry:  as  the 

shear  increased  progressively,  the  magnetic-field  energy  increased 
but  approached  the  energy  of  the  open-field  configuration  asympotically , 
while  in  our  earlier  work  with  the  cylindrical  model,  the  field  energy 
continued  increasing  linearly  past  the  open  field  energy. 

We  spent  some  time  and  energy  with  the  aim  of  determining  whether 
this  effect  is  real  or  whether  it  was  due  to  some  difference  in 
boundary  conditions,  error  in  one  of  the  programs,  etc.  In  order  to 
help  resolve  these  questions,  we  decided  to  simulate  the  cylindrical 


model  by  means  of  the  spherical  code.  This  could  be  done  by  considering 
a ring  dipole  very  close  to  one  of  the  poles  of  the  sphere.  This  model 
substantially  reproduced  the  results  of  the  earlier  cylindrical  calcu- 
lations . 

This  test  indicated  that  there  are  qualitative  differences  between 

the  cylindrical  and  spherical  models.  Figures  3*5  ~ 3.8  show  results 

from  these  runs.  The  R and  Z boundaries  are  conducting  in 

max  max 

Figures  3*5  - 3.6;  the  shear  angles  are  tt  and  2tt  respectively.  In 
Figure  3*5,  the  outermost  field  lines  are  "flattened"  against  the 
boundaries,  particularly  the  R ^ boundary.  This  effect  is  even  more 
pronounced  in  Figure  3-6,  where  the  fields  are  strongly  influenced  by 
the  outer  boundaries. 

Figure  3.7  shows  the  shear  2rr  case  with  the  two  outer  boundaries 
opened  up.  The  run  associated  with  this  figure  was  slow  to  converge, 
requiring  a small  relaxation  parameter  for  stability.  It  is  possible 
that  further  iteration  would  see  the  lines  open  up  further,  although 
we  feel  the  configuration  shown  is  near  the  final  value.  Note  that 
the  fields  open  outwards  only  about  ten  times  the  radius  of  the  neutral 
shear  radius.  This  is  in  marked  contrast  to  the  spherical  case,  where 
they  extended  out  several  hundred  times  the  scale  of  the  neutral  shear 
line  radius.  The  energy  stored  in  the  fields  fell  somewhat  when  the 
boundaries  were  opened  for  the  shear  of  Figure  3.7>  it  was  greater  than 
the  completely  open  field  energy^.  We  cannot  say  whether  this  is  because 
an  appreciable  amount  of  energy  was  carried  across  the  boundaries  when 
the  lines  opened  up,  or  because  the  more  open  state  is  simply  a state  of 
significantly  lower  energy.  Certainly,  the  energy  versus  shear  curve 
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with  conducting  boundaries  increases  almost  linearly  past  the  completely 
open  field  value,  in  agreement  with  our  earlier  result  (1972).  Figure  3.8 
shows  the  field  configuration  when  the  lines  are  completely  open. 

This  simulation  of  the  earlier  cylindrical  code  with  the  spherical 
model  indicates  that  there  is  a physical  difference  between  the  two 
results  and  that  there  are  apparently  two  types  of  force-free  field 
configurations : 

Type  1.  For  a given  magnetic  flux  distribution  on  the  bounding  surface, 
as  the  shear  increases,  the  stored  energy  increases  but  approaches,  as 
a limiting  value,  the  energy  of  the  corresponding  open-field  configuration. 
Type  2.  For  given  magnetic  flux  distribution  on  the  bounding  surface, 
the  stored  energy  increases  indefinitely  as  the  amount  of  shear  increases. 
In  particular,  after  a finite  amount  of  shear,  the  stored  energy  exceeds 
the  energy  of  the  corresponding  open-field  configuration. 

Nevertheless,  it  must  be  emphasized  that  we  have  only  a presumptive 
case  for  this  classification:  it  will  be  essential  to  investigate 

a wider  range  of  models  and,  even  more  important,  to  pin  down  the 
precise  effect  of  the  outer  boundary  condition,  before  this  classifi- 
cation can  be  considered  to  be  established.  It  remains  at  this  time 
a possibility  that  the  apparent  classification  is  an  artifact  of 
subtle  differences  in  boundary  conditions  used  for  different  models. 

The  models  of  cylindrical  topology  which  we  have  considered 
belong  to  Type  2.  The  early  models  of  spherical  geometry,  which  we 
have  calculated  as  part  of  this  program,  belong  to  Type  1.  However, 
when  we  used  the  spherical  geometry  but  compressed  the  field  configu- 
ration to  a small  region  around  one  of  the  poles,  the  behavior  corre- 
sponded to  that  of  Type  2. 


11 


IV.  DISCUSSION 


We  have  found  that  it  is  possible  to  compute  force-free  field 
configurations  without  making  the  restrictive  assumption  (adopted  in 
most  other  work)  that  * = constant.  However,  the  convergence  can 
become  very  slow  if  the  magnetic  field  is  highly  stressed,  and  it  is 
clearly  desirable  to  try  to  find  some  procedures  to  speed  up  convergence. 
It  also  would  be  highly  advantageous  to  pursue  these  calculations  with 
a more  powerful  computer. 

The  most  important  result  of  this  program  is  the  indication  that 
force-free  magnetic-field  configurations  may  be  divided  into  two 
types,  defined  at  the  end  of  Section  III.  The  evolution  of  force-free 
fields  of  Type  1 would  probably  not  lead  to  a solar  flare.  As  the 
stress  increases  progressively,  the  field  expands  and  approaches 
asympotically  the  open  state.  Since  this  open  configuration  builds 
up  slowly,  it  seems  unlikely  that  major  reconnection  would  occur.  It 
seems  more  likely  that  the  field  topology  would  develop  into  a geometry 
with  a thick  current  sheet  which  is  fairly  stable,  such  as  exists  in  a 
coronal  streamer.  However,  it  is  always  possible  that  a certain 
disturbance  (such  as  a shock  wave)  might  trigger  an  instability  even 
in  such  a current  sheet. 

By  contrast,  force-free  fields  of  Type  2 are  highly  important  for 
the  flare  problem.  They  indicate  the  following  sequence  of  events  in 
an  active  region  leading  to  a flare:  a subphotospheric  twisted  flux 

tube  erupts  to  give  a bipolar  field  configuration  in  an  active  region. 
Since  it  is  energetically  favorable  for  the  twist  to  propagate  from 
below  the  photosphere  (where  the  field  is  strong)  to  above  the  photosphere 
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(where  it  is  weaker),  this  will  occur,  leading  to  progressive  twisting 
of  the  force-free  field  configuration  in  the  atmosphere  above  the 
active  region.  When  the  field  is  sufficiently  stressed,  it  contains 
more  energy  than  does  the  corresponding  open-field  configuration.  It 
is  therefore  either  unstable  or  metastable.  The  fact  that  the 
relaxation  calculation  is  stable  (although  only  marginally  so)  indicates 
that  such  a field  configuration  is  physically  stable  against  small 
perturbations.  It  must  therefore  be  metastable,  that  is,  unstable 
against  finite  perturbations.  A finite  perturbation  would  therefore 
trigger  an  eruptive  instability  converting  the  closed  force-free 
configuration  into  the  corresponding  open-field  configuration  with  a 
current  sheet.  The  magnetic-field  reconnection  can  occur  immediately 
in  this  current  sheet  leading  to  a solar  flare. 

To  sum  up,  it  seems  that  the  evolution  of  a force-free  field  of 
Type  1 would  probably  not  lead  to  a flare,  whereas  the  evolution  of 
force-free  fields  of  Type  2 would  almost  certainly  lead  to  a flare. 

It  is  important  to  pursue  the  indication  of  the  existence  of  two 
types  of  force-free  field.  A wider  range  of  models  should  be  calculated, 
on  a more  powerful  computer,  either  removing  the  outer  boundary  to  very 
large  distances,  or  finding  a method  which  does  not  require  an  outer 
boundary.  For  instance,  it  might  be  possible  to  replace  the  radial 
variable  by  1/r.  However,  in  this  case,  it  will  be  essential  to  study 
asymptotic  forms  of  force-free  magnetic-field  configurations  to  insure 
that  we  can  properly  handle  the  behavior  of  the  magnetic  field  in  the 
neighborhood  of  the  singularity  1/r  = 0. 
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FIGURE  CAPTIONS 


Figure  3.1.  Field  configuration  for  a spherical  case  with  a conducting 
equator.  The  neutral  shear  point  occurs  at  about  65°  N. 
latitude.  The  radial  coordinate  is  logarithmic  and  extends 
to  l80  surface  radii.  Maximum  shear  for  this  case  is  2.6 
radians . 


Figure  3*2.  The  same  as  Figure  3.1  with  maximum  shear  tt  radians. 

Figure  3*3*  The  same  as  Figure  3*2  with  maximum  shear  2tt  radians. 

Figure  3«^»  Curve  of  stored  energy  vs  maximum  shear  angle  for  the  case 

of  Figures  3.1  - 3*3.  Shown  also  is  the  approximate  value 
of  the  stored  energy  of  the  completely  open  field  config- 
uration . 


Figure  3«5»  Field  configuration  for  a case  in  which  the  computational 
region  is  restricted  to  be  near  the  pole,  simulating  the 
earlier  cylindrical  model.  The  system  is  nearly  "square" 
(that  is  a cylinder  with  radius  equal  to  length)  with  the 
neutral  shear  point  occurring  at  about  1/7  of  the  outer 
radius.  Maximum  shear  for  this  case  is  tt  radians. 


Figure  3.6.  The  same  as  Figure  3.5  with  maximum  shear  2n  radians. 


Figure  3*7.  The  same  as  Figure  3.6  with  the  Z and  R boundaries 
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opened  to  allow  field  penetration. 


Figure  3*8.  The  same  configuration  as  Figure  3*5  with  completely  open 
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field  lines.  Here  the  pattern  is  current  free  except  along 
the  contour  beginning  at  the  neutral  shear  point,  where  a 
current  sheet  exists. 
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J*  4.*  V-  Ml  ' Vx 


MAXIMUM  SHEAR  AN 


